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Abstract: Small but macroscopic particles - chondrules, higher temperature min- 
eral inclusions, metal grains, and their like - dominate the fabric of primitive 
meteorites. The properties of these constituents, and their relationship to the 
fine dust grains which surround them, suggest that they led an extended exis- 
tence in a gaseous protoplanetary nebula prior to their incorporation into their 
parent primitive bodies. In this paper we explore in some detail the velocities 
acquired by such particles in a turbulent nebula. We treat velocities in inertial 
space (relevant to diffusion), velocities relative to the gas and entrained micro- 
scopic dust (relevant to accretion of dust rims), and velocities relative to each 
other (relevant to collisions). We extend previous work by presenting explicit, 
closed-form solutions for the magnitude and size dependence of these velocities 
in this important particle size regime, and compare these expressions with new 
numerical calculations. The magnitude and size x dependence of these velocities 
have immediate applications to chondrule and CAI rimming by fine dust, and 
to their diffusion in the nebula, which we explore separately. 


1 Background 

The fabric of the most primitive meteorites undoubtedly contains many clues as to their 
origin. While most chondrites are samples of surfaces that have been well worked over by 
impacts and stirring (“regolith breccias”), the dominance of chondrules and like-sized objects 
remains clear. How it came about that most chondrite parent bodies are so dominated by 
particles with such a well-defined range of physical, chemical, and petrographic properties 
remains one of the big puzzles of meteoritics. Since there are relatively few examples of 
anything larger than 0.1-10 mm size particles in most primitive planetesimals, the way such 
particles interact with the gaseous nebula is of prime importance. 

Fe-Mg-Si-0 mineral chondrules, which solidified from a melt, constitute 30-80% of primi- 
tive meteorites. There are a number of extant hypotheses for the formation of the chondrules. 
Most workers in the field believe that chondrules are formed by either localized or nebula 
scale energetic events operating on freely floating precursors of comparable mass, at some 
location or locations in the protoplanetary nebula. However, some still maintain they are 
made in or on primitive bodies, or in collisions between them. In a hybrid scenario, some 
suggest they are formed in shock waves generated by already-formed planetesimals, and thus 
that they are a secondary phenomenon to primary accretion of planetesimals. See eg. Gross- 
man (1989), Grossman et al. (1989), Boss (1996), Connolly and Love (1999), and Jones et 
al. (2000) for reviews of hypotheses on this long-controversial and perennially fascinating 
subject. 
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Another meteorite constituent of great interest are the mineral grains called Ca-Al-rich 
refractory inclusions (CAIs) - so called because their constituent minerals condense out of 
nebula gas at a much higher temperature than do chondrules. These objects are widely 
believed to be direct nebula condensates, and have a complex subsequent thermal history 
which has some similarities to that of chondrules and some differences. There is some in- 
dication from radioisotope ages that CAIs might be ~ 10 6 years older than the chondrules, 
but this remains slightly controversial. They make up 1-10% of primitive meteorites de- 
pending on type, and their size distribution is broader than that of the chondrules. How 
these high-temperature minerals find themselves intimately mixed with lower-temperature 
minerals remains a puzzle. 

It remains unresolved at this time whether the nebula gas was turbulent or laminar dur- 
ing the chondrule era. In previous papers, we have suggested that some of the observed 
properties of chondrules themselves - their typical size and size distribution - can be as- 
sociated with, and easily explained by, the effects of weak nebula turbulence (Cuzzi et al 
1996, 2001). Nevertheless, a consistent end-to-end scenario for formation of primitive bodies 
in this environment, and relying on these processes, is not yet in hand. In this paper, we 
focus on the velocity evolution of this specific class of particles in a weakly turbulent nebula 
as a step towards developing a more complete scenario that operates to produce primitive 
bodies in a similar way across a variety of environments. The velocity evolution is critical for 
our understanding of several important aspects of chondrules and chondrites: (a) the radial 
distribution and redistribution or transport of chondrules and/or CAIs, once formed, before 
their accumulation into parent bodies; (b) The presence of fine grained rims on chondrules, 
CAIs, and other coarse particles in primitive chondrites (Metzler and Bischoff 1996, Brearley 
and Jones 1998); and (c) collision rates and velocities between chondrule-sized particles. The 
main goal of this paper is to provide a theoretical framework within which we can better 
understand mm-to-cm-size particle evolution in general. We accomplish this in sections 1 
(analytical theory) and 2 (supporting numerical calculations). In another paper we apply 
these results to diffusion and dust rimming (Cuzzi 2002b). 

1.1 Particle Velocities in Turbulence 

Astrophysical modeling of the basic physics of particle behavior in fluid flows, laminar or 
turbulent, tends to begin and end with the classic papers by Whipple (1973), Adachi et al. 
(1976), Weidenschilling (1977, 1980), and Volk et al (1980, henceforth VJMR; also Volk et 
al 1978), with important recent updates by Markiewicz et al (1991; henceforth MMV). In 
the fluid dynamics literature, however, the study of particle motions in fluid flows has both a 
long history, and a robust ongoing presence. This history is nicely summarized by Meek and 
Jones (1973). More recent work in the fluids literature is noted in various relevant places 
below. VJMR first developed a useful formalism for calculating the dispersion velocities V p 
(relative to inertial space) and collision velocities V pp (relative to each other) of particles in 
a turbulent nebula. They circumvented the thorny problem of “essential nonlinearity” (cf. 
Meek and Jones 1973) by translating clever physical insights into mathematics and adopting 
a velocity autocorrelation function approach, which we discuss in more detail below. While it 
serves an important internal role in their solutions, neither VJMR nor MMV say much about 
the relative velocity between particles and gas, V P9 . Yet, V pg is the determinant quantity for 
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accretion of rims of fine dust grains by small, macroscopic objects (Paque and Cuzzi 1997, 
Cuzzi et al 1998, Morfill et al. 1998). Our goal in this paper is to quantify V p , V pg , and V pp 
for such particles in a way that extends and focusses the formulation of VJMR and MMV, 
and which allows insights to be gained into the history of chondrules and like-sized particles 
in the protoplanetary nebula. 

In this paper, we determine velocities of all three kinds - V p , V pg , and V pp - with emphasis 
on particles having stopping times t s comparable to the overturn time t v of Kolmogorov 
scale eddies. Particles in this size regime have behavior more complex than tiny “dust” 
grains, which are essentially trapped to the gas flow on all scales. In particular, particles 
with t s = t v are subject to “preferential concentration” by large factors in turbulence, and 
based on some of its apparent fingerprints in the meteorite record, we have suggested a link 
between this process, chondrules, and primary accretion. Specifically, we refer to the fact 
that the typical size and the shape of the size distribution of chondrules are readily explained 
by turbulent concentration (Cuzzi et al 1996, 2001). In a parallel paper (Cuzzi 2002b) we 
explore the possibility that the functional form of V pg might reveal another fingerprint of 
turbulent concentration, and that turbulence might help us understand the puzzling mix of 
CAIs and chondrules in the same meteorites. 

Particles are aerodynamically classified by their Stokes number St, the ratio of their 
stopping time t s to the overturn time of some characteristic eddy. We will make use of 
Stokes numbers defined relative to two different eddy overturn timescales: the Stokes number 
relative to the largest, or integral scale eddy time ti\ Sti = t s /ti, and that defined relative to 
the smallest, or Kolmogorov scale eddy time t v : St n = ts/t^. The overturn time of the largest 
scale eddy tL is generally regarded as the local orbit period. Preferentially concentrated 
particles (chondrules, we believe) have St,, = 1 and StL << 1. For these particles, which 
are smaller than the gas molecular mean free path, the stopping time t s = rp s /cp g , where 
r is particle radius, p s is particle material density, c is the nebula sound speed, and p g is 
the nebula gas density (Weidenschilling 1977). That is, t s and thus both Sti and St v are 
linearly proportional to particle radius. 

1.2 Previous work; the autocorrelation function 

We briefly review and simplify the notation of VJMR and MMV. VJMR assumed a fully 
developed inertial range of turbulence with some largest, or integral scale L and zero smallest 
scale. MMV also adopted the Kolmogorov energy spectrum (as shall we) but correctly 
pointed out that turbulence ceases for scales smaller than the Kolmogorov or inner scale 
r/. Especially for small particles in the chondrule-and-CAI size range, MMV point out that 
this has important implications for V p and V pp , and we will show that the implications are 
important for V pg as well. In a Kolmogorov spectrum, an inertial range of turbulent gas 
kinetic energy extends from the largest or integral scale l — L to the smallest or Kolmogorov 
scale l = 7 ]. Following VJMR, we work in the spatial frequency regime, where k(l) = 2-k/I 
and E(k) = E L (k/ki)~ 5//3 for the Kolmogorov spectrum (note our E(k) is a true energy, and 
is half of VJMR’s P(k)). Then v{k) = {2kE{k))V 2 and t(k) = 1 /(kv(k)) = t L (k/k L )~V 3 . 
As did MMV, we assume E(k) — 0 for k > k v (no turbulent energy at scales smaller than 
the Kolmogorov scale). The mean square turbulent (fluctuating) gas velocity is V g \ thus 
the typical turbulent kinetic energy per unit gas mass is V 2 /2, providing the normalization 
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criterion: 


/ E(k)dk = Vg/2 = [*' E L (k/k L )- & / 3 dk = lE L k L . (1) 

The turbulent gas motions induce fluctuating velocities in the particle population, leading 
to diffusion (V p ), mutual collisions ( V pp ), and motion relative to the local gas (V^ 9 ). 

VJMR derive V p formally by a backwards time integration of the instantaneous acceler- 
ation (their equations 5 and 6): 

V p (t) = tj 1 [ ex P (”(* “ t')/t s )V g {t')dt' (2) 

where V g (t') represents the fluctuating gas velocity history along a particle trajectory (for- 
mally unknown at this point). They proceed by approximating V g (t') as an integral over all 
(independently acting) spatial frequencies k with eddy timescales t k , and approximate the 
contributions as coming from two classes of eddies: “class 1” eddies, with overturn times 
long enough ( t k > t s ) that particles are alw r ays in equilibrium within them, and are primarily 
just advected by their (temporally fluctuating) motions, and “class 3” eddies with overturn 
times too short (t k < t s ) for the particle to come to equilibrium as it passes through them. 
Intermediate, or what might be “class 2” eddies are not treated separately, but simply ab- 
sorbed into the classes on either side. Different simplifications are allowed for each class. 
The boundary between eddy classes 1 and 3 is k *, where"!*, = t s . VJMR show that the class 
3 (small, fast) eddies are negligible for velocity components V p and V pg , but dominate the 
contributions to V pp . We will make use of these results below. 

VJMR first obtain the product < V p (t)V p (t) >=< V p > by integrating backwards over 
two separate time histories. They introduce the gas velocity autocorrelation function for gas 
velocities (in their equation 16) R(t,t';k ) = exp (— \t — t'\/t k ). While they don’t make the 
distinction, the autocorrelation function to be used in this way is properly that determined 
along a particle trajectory (Batchelor 1948, Hinze 1975, Squires and Eaton 1990, Elghobashi 
1991), and is thus a function of t s in general. However, for St L « 1, and at this stage of 
our knowledge, this distinction is not significant (Squires 1990). 

Subsequently, MMV suggested a more general, even if ad hoc, functional form for R(t, t'\ k)\ 

R(t, k)=(l + j ( 3 ) 

with n = 0 or 1. They note that the n — 1 case has more plausible physical behavior (zero 
slope) near t = t 1 than the n = 0 (pure exponential) form assumed by VJMR. 

1.2.1 New results regarding the form of the autocorrelation function, and the 
value of n: 

The selection of n = (0, 1) determines the form of the fluid velocity autocorrelation function 
R(t,t';k). Squires (1990) measured this function directly in his direct numerical simulations 
of turbulence, by following fluid motions alorig the trajectories of a number of particles with 
different St L . In figure 1 we compare the results of Squires (1990) with the predictions 
based on the n = 0 and n = 1 expressions of MMV for R(t,t'\k). Note that, since MMV 
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Figure 1: Autocorrelation function for gas velocities along the trajectory of a 5^ = 1 
particle, as computed directly from our simulations (dotted) and from the simulations of 
Squires (1990), and as calculated using the n — 0 and n = 1 models of MMV. Here, r = t — t‘ 
and is normalized by the large eddy turnover time T^. The n = 1 model is clearly the better 
choice. 

r 

express their autocorrelation function as a function of k, it must be integrated over an energy 
spectrum to compare with the numerical results of Squires (1990). Because Squires (1990) 
only calculated a 1-D autocorrelation function (ie., using only one velocity component), we 
integrated the R(t, t'\ k ) of MMV over a ID energy spectrum (essentially, one-third of the 
total E(k )) (see also Squires and Eaton 1991). It is clear from figure 1 that n = 1 is the better 
choice. This has important implications, primarily for V pg and V pp . In section 2, we directly 
compare V p and V pg calculated in full 3D turbulence using the two alternate autocorrelation 
functions, and again reach the same conclusion. 

1.3 Particle random velocities relative to inertial space 

After some algebra, VJMR derive an expression (their equation 18) for the mean square par- 
ticle fluctuating velocity V p , of which we need only the large, slow (class 1) eddy contribution 
since the small eddy contribution is negligible for St v = 1 particles (we will henceforth drop 
the <> notation on V p , V g , V pg , and V pp , and will merely recall that all are statistical expec- 
tation values based on extensive temporal or spatial averaging). Because of our emphasis 
on particles with St v = 1, we also replace the upper limit of VJMR’s class 1 integral ( k *) 
with the Kolmogorov scale k v . This simplification is, in fact, actually fairly good over the 
entire range of St L « 1, precisely because the contribution of eddies on smaller scales than 
k* (the class 3 eddies) is negligible. That is, the upper limit can be extended from k* to k ^ 
in general for mathematical simplicity without incurring significant error. Mathematically, 
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the upper limit could even be extended to infinity {eg., Volk et al. 1980), but the important 
role of the Reynolds number and of the Kolmogorov scale is then lost. Thus, 

v l ~ 2 / #) rr, dk. (4) 

l k T- t s 

Similarly, the generalized MMV expression for V 2 (their equation 6) can be simplified to 


V 2 

p 





( 1 
\1 + ik/ts 



(5) 


for the particle size regime of interest here. As did VJMR, MMV note that the second 
integral of their equation (6) - the class 3 eddy contribution - is negligible for small particles, 
so we retain only the first integral of their equation (6). We again simplify the upper limit 
of integration in the remaining integral for the nominal Sf, to 1 case where k* ss k v >> k L . 
We validate this by comparing our results with those of MMV (section 1.7). 

The result for V 2 was plotted, but not stated explicitly, by VJMR and MMV (figure 1 in 
both papers), and explicitly derived by Cuzzi et al. (1993; Appendix B): V 2 ~ V 2 /{\ + St L ). 
It is simple to see why V 2 to V 2 in the limit St L « 1 and certainly for St v to 1, since 
<< ik in equations (4) or (5) for nearly all k and overwhelmingly all E{k). This limit is 
appropriate fpr chondrule-and-CAI-sized particles even in the presence of their small vertical 
settling velocity - they diffuse nearly as well as a gas molecule, and do not “settle to the 
midplane” in even a very weakly turbulent nebula (Dubrulle et al. 1995, Cuzzi et al. 1996). 
The implications are discussed in section 3. However, V 2 and V 2 are not exactly equal, 
resulting in a small, but very important, relative energy of motion V pg , giving the velocity 
with which particles move through the gas and encounter tiny (micron-sized) dust grains. 


1.4 Particle velocities relative to the gas 

The average relative velocity magnitude between a particle and the turbulent gas is V pg . 
VJMR make use only of the spatial frequency components of this quantity, which they refer to 
35 V-eK^) (their equation 15). Practically speaking, however, a particle will instantaneously 
sense all eddy contributions as one V pg \ we obtain this by merely integrating VJMR equation 
(15) over k. Considering only the part of the expression relevant for St n to 1 (that for 
k* > k L ), neglecting any systematic velocity, and again letting k * to k v » k L , the second 
integral vanishes and we obtain 

v > 2 C m (^rr) dk - <•> 

For this n = 0 case treated by VJMR, it can be easily verified using equations (4) and 
(6) that 

k + v p 2 = 2 r = y g- ( 7 ) 

Jk L 

However, this useful result is true independent of n. It may also be obtained by Fourier 
transform solution of the forcing equations in temporal frequency ( oj ) space, where the 
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energy spectrum of gas velocity fluctuations E g ( u), particle velocity fluctuations E p (ui), and 
relative velocity fluctuations E pg (i <j) are related by 

E p (uj) = E g (u)/{ 1 + t 2 u 2 ) and E pg (u) = t]u) 2 E p (u ) . (8) 


This approach can be traced to Csanady (1963); it is also described by Hinze (1975, chapter 
5), Meek and Jones (1973), and Squires (1990, sections 4.2 and 4.5.1)). The E p solution 
was also derived in this way by Cuzzi et al. (1993, Appendix B). It is also clear then that 
E pg (u>) + E p (u) = E g (uj ), essentially the same result as equation (7) above. Finally, we have 
directly verified equation (7) in our numerical simulations. 

Using this general relationship, we can extend the results of MMV to obtain V pg for their 
more generalized gas velocity autocorrelation functions (they only present results for V 2 ). 
That is, using equations (1) and (5), 


V 2 = V - V 

P9 9 P 


'■L 


ln E(k)(-r~-T+'dk = 2 f k ' E(k) 

ftL i 


n-j-1 


l + t k /t s 


dk. 


(9) 


We will use equations (5) and (9), with assumed inertial range expressions for E(k), to 
derive analytical expressions for V p and V pg of hypothetically “chondrule-like” ( ie St v « 1) 
particles as functions of their size and the turbulent Reynolds number. 


1.5 


Relative velocities between 


particles of similar sizes 


Expressions for V pp (VJMR Appendix C and equation 19; MMV equations 7 and 8) are more 
cumbersome, but respond nicely to certain simplifying assumptions. The full expression for 
V pp for two particles of equal size is (changing notation slightly from MMV equation 9, and 
allowing for a finite Kolmogorov scale): 




t s \ f x , nt s h(x) 

tT+tf) + 


dk, 


( 10 ) 


where g(x) = tan -1 (x)/x and h(x) = 1/(1 + X 2 )- T he parameter x of VJMR and MMV is 
small in our regime of interest: 


Vrel(k)t,{kt k ) = Vrel(k)t a _ Vreijk) 

t s + t k v{k){t s +t k ) ~ 2v(k) 


( 11 ) 


since in the very limited range of k over which the integral is done, t s « t k . 1 In fact x << 1 
over most of the integral where t s << t k) so the functions g(x) and h(x) are ~ 1 or perhaps 
as small as a fraction of order unity; thus 


V 2 

VP 


fkrj 

/ E(k) 

Jk* 


1 - 


t. 


is + tk 


n+11 


dk 


r k i 

= 4 E(k) 

Jk’ 


1 

1 + tk/t s 


Tl+ 1 


dk. (12) 


The integrand is identical to that for V 2 , but the integral has different limits which make it 
clear that only the eddies faster than t s can perturb identical particles into having incoherent 
relative velocities. 

Rn the above equation, the mathematical generalization of V pg by VJMR and MMV to its k — th com- 
ponents Vrei(k) momentarily reappears. However, it is true in general, at any spatial frequency, that the 
particle-gas relative velocity is less than, or at most equal to, the gas velocity itself. 
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1.6 Scaling relations 

Recall that for the gas, 

t k = l(k)/v(k) = {L/V g )(k/k L )~ 2/3 = t L (k/k L )~ 2 ' 3 (13) 

(Cuzzi et al. 2001). In equation (13) we have made the usual identification of V g with the 
largest scale eddy L. For the particles, 

= st L = (k,! k L y 2 / 3 (14) 

and 

~ = ( k/k s ) 2/ 3 = j-{k/k L ) 2/ 3 = St L (k/ k L ) 2/3 , (15) 

Note that if we restrict our attention to particles with St v = t s /t v « 1, then their Stokes 
number referred to the integral scale automatically becomes 

St L = t s /t L = t v /t L = (V*l)- 2/3 = (i?e 3/4 ) _2/3 = Re~ 1/2 . (16) 

The last substitution of (k v /ki) = Re 3 / 4 , where Re = LV g /v is the flow Reynolds number, 
with v being the molecular kinematic viscosity, is a direct consequence of the definitions of 
the Kolmogorov scale, the energy dissipation rate, and tbe Reynolds number (Tennekes and 
Lumley 1972). This relation can be obtained without any reference at all to the Kolmogorov 
spectrum but merely using scaling arguments relating to t L and t v 2 Re is related to astro- 
physical “a” -models of the protoplanetary nebula by Re = acH/u with c = sound speed 
and H = nebula vertical scale height (Cuzzi et al. 2001). 


1.7 Final expressions for V pg and V pp 

Substituting the scaling relations from above for t s /t k , equation (9) for V pg becomes 




St L 

st L + (k/k L )-v* 


n + l 


dk. 


(17) 


We use the normalization (equation 1) to write E{k) = (V 2 /3 k L )(k/k L )~ b / 3 , and change 
integration variable to x = k/k L , leaving 



2V 2 rRe 3/4 / St i \ n+1 

3 Ji Ut L + x- 2 / 3 J 


x-V 3 dx. 


(18) 


where in the upper limit we have substituted k v /k L = Re 3/4 from the scaling relations. 
Closed form solutions for equation (18) can be obtained for n = 0 or 1. For example, for 
n = 1 the result of the integral is 


St L 

1 

— 1/2 

St 2 L {Re 1 ! 2 - 1) 

.1 + Stix 2 / 3 . 

9 

! Re*/ 4 

l(5t L + l)(F^Re 1 / 2 + l)_ 


(19) 


Let the energy dissipation rate be e. Then e — V g /t /, — V g jL where the first expression defines t k 
and the last expression defines L. Also t„ = We) 1/2 and tj = (i^/e) 1 / 4 . Solving gives t L /t„ = Re 1 / 2 and 
i i/L = Re 3/i . 
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For n = 0 the result of the integral is: 



Sti In 


(Re^(l + St L )Y 
V Re l / 2 St L + 1 ) ' 


(20) 


These results make it quite easy to predict both the magnitude and the St v dependence of 
V pg for arbitrary nebula turbulent intensity. 

We solve equation (12) for V pp in a similar fashion to the solution for V pg above, to obtain 

fOr Tl “ l > - . rj In q In 

4V ? rkn/kL f 2St L X~ 7 / 3 + X~ 9 / 3 
Jk ( 


V 2 = 

pp 


3 Jk(t,)/k L \Stl + 2St L x~ 2 / 3 + x" 4 / 3 


dx. 


(21) 


As before, the upper integration limit is k p /k L = Re 3/i . For the lower limit, k*/k L = 
k(t s )/k L - ( t s /t L )~ 31 2 = Stl* /2 from the scaling relations. The closed form analytic solution 
of this integral is: 


V 2 = 2V„ 

PP 9 


' x- 2 / 3 

1 + St L x 2 / 3 


iSi 


- 3/2 


J rte 3 / 4 


= 2V 2 

9 


St L 


1 


St L Re + Re'/ 2 \ 


( 22 ) 


The n = 0 form of the solution is somewhat less useful, and we note it without expanding it 
as it will not be used further. 


V 2 = 2V 2 

r fin r n 


PP 


St L In 


T +St L x 2 >‘ 3 \ 1 


1 iie 3 / 4 


X 2 / 3 


X 2 / 3 


St 


- 3/2 


1.8 Detailed comparisons with the models of Markiewicz et al. 

In addition to developing the analytical expressions discussed and applied in the paper, 
we also developed a detailed numerical model following the prescriptions of MMV exactly 
(but with a generalized turbulent energy spectrum). This was needed both to evaluate 
their theoretical approach in the context of our numerical simulations of turbulence (section 
2), which have a non-Kolmogorov spectrum and low Reynolds number compared to nebula 
applications, and to assess the validity of our analytical approximations. The numerical 
model of MMV is no longer in active use (W. Markiewicz, personal communication 2002), 
so we digitized their V pp results (their figure 5) to facilitate comparisons. As seen in figure 
4, our full numerical model for V pp (solid curves) agrees very well with their results for V pp 
(long dashed curves). In figure 2 we also show our results for V pg , not presented by VJMR 
or MMV, as obtained by integrating MMV equation 4 over all spatial frequencies. Note that 
we, and MMV, both use the appropriate form of R(t,t']k ) (te., that for the correct choice 
of n; section 1.2) for these calculations. 

The most striking feature of the results, first noted by MMV, is that V pp very quickly falls 
to zero for particles with Str, < 1 (he. St L < Re~ l/2 , as shown in the scaling relations of sec- 
tion 1.6 above) because there is no more energy in faster eddies to provide relative velocities 
to such particles. This does not happen to V pg , because eddies on all scales contribute. Also 
note that V p and V pp decrease for large particles ( St L > 1), as fewer eddies can effectively 
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St L = t./L L 


Figure 2: Comparison of our numerical version of the full MMV model for n — 0 (light 
curves) and n = 1 (heavy curves), along with digitized results for V pp from MMV (dashed 
curves, their figure 5, for n = 1). Three different Re are shown: (a) 10 4 , (b) 10 7 , and (c) 
10 9 . The dash-dot curves are for V p , which has the same shape for all three Re. V pg is shown 
in the two sets of dotted curves and V pp in the two sets of solid curves. Note that the n = 0 
values of V pg (light dotted curves) are considerably (3-4 times) higher than the preferred 
n = 1 values (heavy dotted curves), and the St^-dependence of V pg , for n = 0, never gets 
much above 0.5, whereas for n = 1 a linear dependence is seen for Sti < /?.e -1 / 2 . As in 
figures 4 and 5, vertical hash marks correspond to Sti = for the three values of Re. 


10 




couple to particles with such long stopping times. Naturally, V pg simply approaches V g for 
these large particles. 

Upon comparing our original analytical results (equations 19 and 22) with our full nu- 
merical model and the MMV results, we found some small quantitative discrepancies at the 
order unity level, as might be expected. The responsible approximations were easily identi- 
fied. First, we approximated the boundary between class 1 and class 3 eddies by t s = t(k*) 
rather than the more complete equation 9 of VJMR and equation 4 of MMV, which ob- 
tains the relevant eddy frequency in the moving frame of the particle and involves V ret (k). 
Comparison of the two criteria revealed that, to a very good approximation, the criterion 
t s = t(k*) gives a value of k* that is too large by a factor close to 2 (figure 3). So, after 
this “calibration”, we merely decrease the lower limit of integration in our equation (22) by 
a factor of 2. Second, even after this correction, our values of V pp are about 20% high. This 
is easily ascribed to our approximation that g(x) and h(x) are equal to unity throughout 
the entire range of k; in fact, they are tens of percent smaller than unity over some part 
of this range, depending on the value of St L . Empirically, this is corrected by multiplying 
our analytical expression for V pp by a constant factor of 0.8. With these two simple adjust- 
ments, each correcting a known oversimplification, our analytical expression for V pp achieves 
very good good agreement with the MMV results, and with our own full numerical model, 
over the relevant range of St L < 0.1 or so. There appears to be no reason to make such 
refinements to our analytical expression for V pg (equatibn 19), because our approximations 
are better justified and the agreement with MMV acceptable. 


1.9 Numerical refinements to the model 

With insights gained from comparison of our numerical and analytical models, we have made 
two small adjustments to equation (22) for V pp which correct for two of our approximations. 
Equation (22) is multiplied by a factor of 0.8, and the upper integration limit (St L 7 ) is 
divided by two, so the first term in the final expression changes from St L / 2 to St L f 1.03 « 
St L . The approximations entering into our expression for V pg are better, so no correction is 
applied. The final equation for V pp is then 



1.6 V/ 



1 

St L Re + ReV 2 . ' 


(23) 


The results of equations (19) and (23) (the preferred and adjusted n = 1 forms), normalized 
by V g , are shown in figure 4 for the same three values of Re as in MMV, and in closeup 
form in figure 5. 

As shown by MMV (their figure 2), and as seen previously in our figure 2, the falloff 
of V pp is extremely steep for St n < 1 (be. Sti < Re~ 1 ^ 2 as shown in the scaling relations 
of section 1.6 above) because there is no more energy in faster eddies to provide relative 
velocities to such particles. 


1.10 Simplification of analytically determined velocity expressions: 

Equations (19) and (23) - for the preferred n = 1 case - are readily simplified in different 
limits of interest. It is simply shown by retaining leading terms that equation (19) for V pg 
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Figure 3: Correction of our approximation k* w St by a factor of two (dash-dot line) 
which brings it into excellent agreement (in our range of validity St L < 0.1) with the ex- 
act numerical solution for k * , shown for Re = 10 4 , 10 7 , and 10 9 , computed using the full 
VJMR/MMV expression. Only very close to St v = 1 does our approximation deviate slightly; 
notice the tiny tail at St L = 6 x 10~ 3 , k* = 10 3 , which is the Kolmogorov scale for Re = 10 4 . 
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Figure 4: V pg (St L ) (dotted; equation 19) and V pp (St L ) (solid; equation 23) for Re = 10 4 (a), 
10 7 (b), and 10 9 (c). The digitized results of MMV (theif figure 5) for V pp , for the same three 
values of Re, are shown by the dashed lines. Our V pp expression is invalid for St L > 0.1 or 
so (see text). 



st L 


Figure 5: A closeup plot of V pg (dotted), V pp (solid), and the digitized MMV results (their 
figure 5) for V pp (dashed, see Appendix) all for the Re = 10 7 case. The dash-dot line has 
slope 1/2. The vertical short dashed line indicates Str, = 1, where St L = Re~ l/2 \ here, 
a St° L 75 . 
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results in three separate regimes: V pg « V g for St L > 1, oc 5t[ /2 for .Re -1 / 2 < « 1, 

and VJ, 9 oc StiRe 1 ^ for Sti < Re~ 1 ^ 2 . This is confirmed by inspection of figures 2 and 4. In 
the special case of St v = 1, or Sti = Re -1 / 2 , equation (19) reduces directly to 

/?g — 1/4 / iy \ 1/4 

V' s (S(, = 1) = = ca'/< (— ) , (24) 

where we have substituted V g = ca 1/2 (Cuzzi et al. 2001). This Re-dependence, which also 
applies for St v < 1 in general, quite naturally explains a result we obtained empirically from 
our numerical models over a range of Re much smaller than nebula values, namely that 
Vpg/Vg c* R e 1//4 (Cuzzi et al. 1998). By contrast, it is similarly shown from equation (20) 

that the St L -dependence of V pg for the older n = 0 case continues the St / 2 dependence to 
arbitrarily small St L . 

These results are also consistent with arguments in Cuzzi et al. (1993, Appendix B; 
A. Dobrovolskis, personal communication). Expand and time-average the instantaneous 
quantity < (V p - V g ) 2 > to obtain < V pg V pg >=< V p V p > + < V g V g > -2 < V p V g >. 
Substituting from Cuzzi et al. (1993, equation Bll) we find < V p V p >=< V p V g >=< V g V g > 

/(I + St L ), leading to V pg = (St L /{ 1 + St L )) l/ 2 V g , which reaches the same limits as equation 
(19) except for particles with t s < t v , or < 1, because the integral in its derivation 
(equation Bll of Cuzzi et al. 1993) extends to infinite eddy frequency. 

Thus, unless t s < tjj ( Sti, < Re 1//2 ), the particle-gas relative velocity in turbulence is 
generally proportional to ^/Stl for small Sh- The steeper dependence of V pg on St L and 
St v is restricted (in turbulence) to particles with St g < 1. That is, evidence for a more 
nearly linear dependence of V pg on r , if the environment was turbulent, would imply that 
the particles in question were St v < 1 particles. This neW result derives directly from the 
use of the n = 1 gas velocity autocorrelation function. The primary qualitative change is in 
the particle size dependence of V pg for particles with St^ < 1. We address the significance 
of this in more detail in a forthcoming paper (Cuzzi 2002b). 

Finally, using equation (23) for V pp , we get 

V^Str, = 1) = VoIVgRe- 1 ' 4 = 1.26V pg , (25) 

where we used equation (24) for V pg . 

2 Comparison with numerical results 

In this section we compare numerical results from our full 3D Lagrangian particle-gas model 
(Hogan et al. 1999) with full numerical calculations using our implementation of MMV 
(sections 1.8-1. 9). We present particle velocities relative to the computational box (V p ), and 
relative to the local fluid velocity (V^ ff ), as obtained from our simulations. These velocities 
are defined as RMS spatial averages over all particles in a single snapshot, or V = (< ( V x — < 
h: >) 2 > + < (h,— < V y >) 2 > + < ( V 2 — < V z >) 2 >) ] / 2 , where V represents V p or V pg at 
the location of each particle, and <> is the averaging operator < ... >= Sil r ’i(- -)/ J %>’ where 
N p is the number of particles in a single snapshot. Of course, < V > is very close to zero 
for both these quantities since there is no mean flow in our simulations. 
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- MMV: n = 0 
• MMV: n= 1 


0 0 1 0.2 0.3 

st L 

Figure 6: V p vs. St L obtained from our direct simulations compared with MMV predictions 
for models n = 0 (solid line) and n = 1 (dashed line). AH velocities have been normalized by 
the RMS fluid velocity V g . Results are shown for three different Re; the St L values for each 
point are defined relative to a large eddy time based on energy dissipation 2 , which varies 
with Re for our numerical calculations. When they are defined relative to a constant large 
eddy time, as are our analytical models and the MMV models shown in figure 2, points 
and models for all Re collapse onto the same curve as seen in figure 2. The n = 1 MMV 
prediction is clearly a better fit to the numerically simulated velocities, regardless of the 
choice of normalization timescale. 

This spatial averaging approach is equivalent to the temporal averaging implicit in the 
MMV model, because of the ergodic principle that equates temporal and spatial averaging 
under suitable conditions. In our case, the conditions are satisified because our integral 
length scale L is small compared to the spatial period of the computational domain, for all 
Re. 

The case of V pp is more complicated, as the results depend on the proximity region chosen 
for “neighboring” particles. For the most useful comparisons with the predictions of MMV 
and VJMR, and with the expected uses of this quantity in mind, the region over which 
particle neighbors are selected should be as small as possible - less than rj certainly - and 
here we run into sampling errors. Perhaps most important, the deviation of our model energy 
spectrum from a Kolmogorov spectrum is significant (eg. Squires and Eaton 1990), and V pp is 
much more sensitive to the details of the high-spatial-frequency end of the energy spectrum 
than either V p or V pg . Since the main purpose of these calculations is to verify numerically 
the preference for the n = 1 autocorrelation function in an independent way from the direct 
comparison shown in figure 1, and because this case is already well made by the V v and V pg 
plots, we present no comparisons for V pp . 

Figures 6-9 show that the n = 1 autocorrelation function provides a much better 
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Figure 7: V p , as in figure 6, but plotted against St v . 



st L 


Figure 8: V pg obtained from our direct simulations, compared with predictions of the MMV 
models with n = 0 (solid line) and n = 1 (dashed line) vs St All curves have been 
normalized by the RMS fluid velocity V g : 





Figure 9: The same data shown in figure 8, but plotted against St 

fit to both V p and V pg than the n = 0 version. For V pg , the Tits of the MMV theory to 
our simulations are less perfect than for V p . We can see several possible explanations for 
this. For instance, the mathematically simple form adopted for the n = 1 autocorrelation 
function is not a perfect fit to the actual numerically determined one (figure 1), by about the 
correct fractional amount. Also, we have emphasized that the correct velocity autocorrelation 
function to use is that along a particle trajectory (Meek and Jones 1973), and this function 
is actually somewhat size dependent even over the range St v ~ 1 (see, eg., Squires 1990, 
figure 4-23). Finally, because of the deviation of our turbulent kinetic energy spectrum from 
an inertial range, some of the definitions of eddy times used in the MMV theory might 
be inappropriate. It would not be surprising for V pg to be more sensitive to these small 
deviations than V p (compare figures 6 and 8, or 7 and 9). In spite of the small deviations in 
V pg , the combination of the direct comparisons of the autocorrelation functions themselves 
(figure 1), and the comparison of the velocities derived using them (figures 6-9) makes it 
clear that the n = 1 autocorrelation function is the best choice. 


3 Summary and conclusions 

We present theoretical and numerical results which describe the turbulence-driven velocities 
of particles in the St L « 1 size regime which might characterize chondrules and similar 
sized particles. We numerically verify the general approach of VJMR as modified by MMV, 
and verify in two different ways the intuitive preference of MMV for an n = 1 gas velocity 
autocorrelation function - at least along the trajectories of St v « 1 particles. We find 
theoretically that the n = 1 autocorrelation function leads to a particle-gas relative velocity 
function that approaches linear dependence on particle size for particles in the St^ « 1 
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regime, and becomes and remains linear for arbitrarily small sizes. This is quite a different 
result than predicted by the original VJMR n = 0 expressions. We derive simple analytic 
expressions for V p , V pg , and V pp (the latter, for comparable size particles only) for arbitrary 
levels of nebula intensity, as characterized by its Reynolds number Re or its corresponding 
“a”. In a separate paper (Cuzzi 2002b) we will present some implications of these results 
for meteoritics. 
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Table I. List of symbols 


Parameter 

c 

C 

m 

H 

k 

kfj 

L 

r 

R 

Re 

St L 

st„ 

ts 

tk 

tL 

tr] 

y s 

Vv 

v P9 

Vpp 

a 

t 

V 

v 

Vrp 

U) 

Pg 

Ps 


Definition 

gas molecule thermal speed 

particle concentration factor 

turbulent kinetic energy at wavenumber k 

nebula vertical scale height 

eddy wavenumber 

wavenumber of largest eddy 

wavenumber of Kolmogorov scale eddy 

integral or largest scale in turbulent energy spectrum 

particle radius 

gas velocity autocorrelation function 

flow Reynolds number 

Stokes number relative to largest eddy 

Stokes number relative to Kolmogorov scale eddy 

stopping time of particle due to gas drag 

overturn time of eddy with wavenumber k 

overturn time of largest eddy 

overturn time of Kolmogorov scale eddy 

gas turbulent velocity (large eddy) 

particle random velocity in inertial space 

relative velocity between particles and gas 

relative velocity between particles r 

nebula viscosity parameter; Re = acH/u 

dissipation of turbulent kinetic energy 

Kolmogorov scale 

molecular kinematic viscosity 

turbulent kinematic viscosity 

eddy temporal frequency 

gas mass density 

solid particle mass density 
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Abstract: Small but macroscopic particles - chondrules, higher temperature min- 
eral inclusions, metal grains, and their like - dominate the fabric of primitive 
meteorites. The properties of these constituents, and their relationship to the 
fine dust grains which surround them, suggest that they led an extended exis- 
tence in a gaseous protoplanetary nebula prior to their incorporation into their 
parent primitive bodies. In this paper we explore in some detail the velocities 
acquired by such particles in a turbulent nebula. We treat velocities in inertial 
space (relevant to diffusion), velocities relative to the gas and entrained micro- 
scopic dust (relevant to accretion of dust rims), and velocities relative to each 
other (relevant to collisions). We extend previous work by presenting explicit, 
closed-form solutions for the magnitude and size dependence of these velocities 
in this important particle size regime, and compare these expressions with new 
numerical calculations. The magnitude and sizer* dependence of these velocities 
have immediate applications to chondrule and CAI rimming by fine dust, and 
to their diffusion in the nebula, which we explore separately. 

1 Background 

The fabric of the most primitive meteorites undoubtedly contains many clues as to their 
origin. While most chondrites are samples of surfaces that have been well worked over by 
impacts and stirring (“regolith breccias”), the dominance of chondrules and like-sized objects 
remains clear. How it came about that most chondrite parent bodies are so dominated by 
particles with such a well-defined range of physical, chemical, and petrographic properties 
remains one of the big puzzles of meteoritics. Since there are relatively few examples of 
anything larger than 0.1-10 mm size particles in most primitive planetesimals, the way such 
particles interact with the gaseous nebula is of prime importance. 

Fe-Mg-Si-0 mineral chondrules, which solidified from a melt, constitute 30-80% of primi- 
tive meteorites. There are a number of extant hypotheses for the formation of the chondrules. 
Most workers in the field believe that chondrules are formed by either localized or nebula 
scale energetic events operating on freely floating precursors of comparable mass, at some 
location or locations in the protoplanetary nebula. However, some still maintain they are 
made in or on primitive bodies, or in collisions between them. In a hybrid scenario, some 
suggest they are formed in shock waves generated by already-formed planetesimals, and thus 
that they are a secondary phenomenon to primary accretion of planetesimals. See eg. Gross- 
man (1989), Grossman et al. (1989), Boss (1996), Connolly and Love (1999), and Jones et 
al. (2000) for reviews of hypotheses on this long-controversial and perennially fascinating 
subject. 
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Another meteorite constituent of great interest are the mineral grains called Ca-Al-rich 
refractory inclusions (CAIs) - so called because their constituent minerals condense out of 
nebula gas at a much higher temperature than do chondrules. These objects are widely 
believed to be direct nebula condensates, and have a complex subsequent thermal history 
which has some similarities to that of chondrules and some differences. There is some in- 
dication from radioisotope ages that CAIs might be ~ 10 6 years older than the chondrules, 
but this remains slightly controversial. They make up 1-10% of primitive meteorites de- 
pending on type, and their size distribution is broader than that of the chondrules. How 
these high-temperature minerals find themselves intimately mixed with lower-temperature 
minerals remains a puzzle. 

It remains unresolved at this time whether the nebula gas was turbulent or laminar dur- 
ing the chondrule era. In previous papers, we have suggested that some of the observed 
properties of chondrules themselves - their typical size and size distribution - can be as- 
sociated with, and easily explained by, the effects of weak nebula turbulence (Cuzzi et al 
1996, 2001). Nevertheless, a consistent end-to-end scenario for formation of primitive bodies 
in this environment, and relying on these processes, is not yet in hand. In this paper, we 
focus on the velocity evolution of this specific class of particles in a weakly turbulent nebula 
as a step towards developing a more complete scenario that operates to produce primitive 
bodies in a similar way across a variety of environments. The velocity evolution is critical for 
our understanding of several important aspects of chondrules and chondrites: (a) the radial 
distribution and redistribution or transport of chondrules and/or CAIs, once formed, before 
their accumulation into parent bodies; (b) The presence of fine grained rims on chondrules, 
CAIs, and other coarse particles in primitive chondrites (Metzler and Bischoff 1996, Brearley 
and Jones 1998); and (c) collision rates and velocities between chondrule-sized particles. The 
main goal of this paper is to provide a theoretical framework within which we can better 
understand mm-to-cm-size particle evolution in general. We accomplish this in sections 1 
(analytical theory) and 2 (supporting numerical calculations). In another paper we apply 
these results to diffusion and dust rimming (Cuzzi 2002b). 

1.1 Particle Velocities in Turbulence 

Astrophysical modeling of the basic physics of particle behavior in fluid flows, laminar or 
turbulent, tends to begin and end with the classic papers by Whipple (1973), Adachi et al. 
(1976), Weidenschilling (1977, 1980), and Volk et al. (1980, henceforth VJMR; also Volk et 
al. 1978), with important recent updates by Markiewicz et al. (1991; henceforth MMV). In 
the fluid dynamics literature, however, the study of particle motions in fluid flows has both a 
long history, and a robust ongoing presence. This history is nicely summarized by Meek and 
Jones (1973). More recent work in the fluids literature is noted in various relevant places 
below. VJMR first developed a useful formalism for calculating the dispersion velocities V p 
(relative to inertial space) and collision velocities V pp (relative to each other) of particles in 
a turbulent nebula. They circumvented the thorny problem of “essential nonlinearity” (c/. 
Meek and Jones 1973) by translating clever physical insights into mathematics and adopting 
a velocity autocorrelation function approach, which we discuss in more detail below. While it 
serves an important internal role in their solutions, neither VJMR nor MMV say much about 
the relative velocity between particles and gas, V pg . Yet, V pg is the determinant quantity for 
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accretion of rims of fine dust grains by small, macroscopic objects (Paque and Cuzzi 1997, 
Cuzzi et al 1998, Morfill et al. 1998). Our goal in this paper is to quantify V p , V pg , and V pp 
for such particles in a way that extends and focusses the formulation of VJMR and MMV, 
and which allows insights to be gained into the history of chondrules and like-sized particles 
in the protoplanetary nebula. 

In this paper, we determine velocities of all three kinds - V p , V pg , and V pp - with emphasis 
on particles having stopping times t s comparable to the overturn time t n of Kolmogorov 
scale eddies. Particles in this size regime have behavior more complex than tiny “dust” 
grains, which are essentially trapped to the gas flow on all scales. In particular, particles 
with t s = t v are subject to “preferential concentration” by large factors in turbulence, and 
based on some of its apparent fingerprints in the meteorite record, we have suggested a link 
between this process, chondrules, and primary accretion. Specifically, we refer to the fact 
that the typical size and the shape of the size distribution of chondrules are readily explained 
by turbulent concentration (Cuzzi et al 1996, 2001). In a parallel paper (Cuzzi 2002b) we 
explore the possibility that the functional form of V pg might reveal another fingerprint of 
turbulent concentration, and that turbulence might help us understand the puzzling mix of 
CAIs and chondrules in the same meteorites. 

Particles are aerodynamically classified by their Stokes number St, the ratio of their 
stopping time t s to the overturn time of some characteristic eddy. We will make use of 
Stokes numbers defined relative to two different eddy overturn timescales: the Stokes number 
relative to the largest, or integral scale eddy time t^: Sti = t s jti , and that defined relative to 
the smallest, or Kolmogorov scale eddy time t v \ St v = t s /t^- The overturn time of the largest 
scale eddy is generally regarded as the local orbit period. Preferentially concentrated 
particles (chondrules, we believe) have St ^ = 1 and St l << 1. For these particles, which 
are smaller than the gas molecular mean free path, the stopping time t s — rp s /cp g , where 
r is particle radius, p s is particle material density, c is the nebula sound speed, and p g is 
the nebula gas density (Weidenschilling 1977). That is, t s and thus both St L and St v are 
linearly proportional to particle radius. 

1.2 Previous work; the autocorrelation function 

We briefly review and simplify the notation of VJMR and MMV. VJMR assumed a fully 
developed inertial range of turbulence with some largest, or integral scale L and zero smallest 
scale. MMV also adopted the Kolmogorov energy spectrum (as shall we) but correctly 
pointed out that turbulence ceases for scales smaller than the Kolmogorov or inner scale 
rj. Especially for small particles in the chondrule-and-CAI size range, MMV point out that 
this has important implications for V p and V pp , and we will show that the implications are 
important for V pg as well. In a Kolmogorov spectrum, an inertial range of turbulent gas 
kinetic energy extends from the largest or integral scale / = L to the smallest or Kolmogorov 
scale / = r). Following VJMR, we work in the spatial frequency regime, where k(l) — 2 tt/1 
and E(k ) = Ei(k/ki )~ 5 / 3 for the Kolmogorov spectrum (note our E(k) is a true energy, and 
is half of VJMR’s P{k)). Then v(k) — (2 kE(k)) 1 / 2 and t{k ) = l/(kv(k)) = t L {k/k L )~ 2 P . 
As did MMV, we assume E(k ) = 0 for k > k p (no turbulent energy at scales smaller than 
the Kolmogorov scale). The mean square turbulent (fluctuating) gas velocity is V g \ thus 
the typical turbulent kinetic energy per unit gas mass is V 2 /2, providing the normalization 
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criterion: 


P E(k)dk = lf/2 = P E L (k/k L )-V 3 dk = \E L k L . (1) 

Jki Z 

The turbulent gas motions induce fluctuating velocities in the particle population, leading 
to diffusion (VJ,), mutual collisions {V pp ), and motion relative to the local gas (Vj, ff ). 

VJMR derive V v formally by a backwards time integration of the instantaneous acceler- 
ation (their equations 5 and 6): 

V P {t) = tj 1 f exp (—(t — t')/t s )V g {t')dt' (2) 

J o 

where V g (t') represents the fluctuating gas velocity history along a particle trajectory (for- 
mally unknown at this point). They proceed by approximating V g (t') as an integral over all 
(independently acting) spatial frequencies k with eddy timescales t k , and approximate the 
contributions as coming from two classes of eddies: “class 1” eddies, with overturn times 
long enough (t k > t s ) that particles are always in equilibrium within them, and are primarily 
just advected by their (temporally fluctuating) motions, and “class 3” eddies with overturn 
times too short (t k < t s ) for the particle to come to equilibrium as it passes through them. 
Intermediate, or what might be “class 2” eddies are not treated separately, but simply ab- 
sorbed into the classes on either side. Different simplifications are allowed for each class. 
The boundary between eddy classes 1 and 3 is k*, wherel** = t s . VJMR show that the class 
3 (small, fast) eddies are negligible for velocity components V p and V pg , but dominate the 
contributions to V pp . We will make use of these results below. 

VJMR first obtain the product < V p {t)V p {t) >=< > by integrating backwards over 

two separate time histories. They introduce the gas velocity autocorrelation function for gas 
velocities (in their equation 16) R(t,t'] k ) = exp(-|f - t'\/t k ). While they don’t make the 
distinction, the autocorrelation function to be used in this way is properly that determined 
along a particle trajectory (Batchelor 1948, Hinze 1975, Squires and Eaton 1990, Elghobashi 
1991), and is thus a function of t. s in general. However, for St L « 1, and at this stage of 
our knowledge, this distinction is not significant (Squires 1990). 

Subsequently, MMV suggested a more general, even if ad hoc , functional form for R(t, t ' ; k ): 

Jl(t,t';k)=f 1 + (3) 

with n = 0 or 1. They note that the n = 1 case has more plausible physical behavior (zero 
slope) near t = t' than the n = 0 (pure exponential) form assumed by VJMR. 

1.2.1 New results regarding the form of the autocorrelation function, and the 
value of n: 

The selection of n = (0, 1) determines the form of the fluid velocity autocorrelation function 
R.(t,t';k). Squires (1990) measured this function directly in his direct numerical simulations 
of turbulence, by following fluid motions along the trajectories of a number of particles with 
different St L . In figure 1 we compare the results of Squires (1990) with the predictions 
based on the n = 0 and n = 1 expressions of MMV for R(t,t'\k). Note that, since MMV 
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Figure 1: Autocorrelation function for gas velocities along the trajectory of a St v = 1 
particle, as computed directly from our simulations (dotted) and from the simulations of 
Squires (1990), and as calculated using the n = 0 and n = 1 models of MMV. Here, r = t-t' 
and is normalized by the large eddy turnover time T^. The n — 1 model is clearly the better 
choice. 

r 

express their autocorrelation function as a function of k, it must be integrated over an energy 
spectrum to compare with the numerical results of Squires (1990). Because Squires (1990) 
only calculated a 1-D autocorrelation function ( ie using only one velocity component), we 
integrated the R(t,t'\k) of MMV over a ID energy spectrum (essentially, one-third of the 
total E(k )) (see also Squires and Eaton 1991). It is clear from figure 1 that n = 1 is the better 
choice. This has important implications, primarily for V pg and V pp . In section 2, we directly 
compare V p and V pg calculated in full 3D turbulence using the two alternate autocorrelation 
functions, and again reach the same conclusion. 

1.3 Particle random velocities relative to inertial space 

After some algebra, VJMR derive an expression (their equation 18) for the mean square par- 
ticle fluctuating velocity V p , of which we need only the large, slow (class 1) eddy contribution 
since the small eddy contribution is negligible for St^ — 1 particles (we will henceforth drop 
the <> notation on V p , V g , V pg , and V pp , and will merely recall that all are statistical expec- 
tation values based on extensive temporal or spatial averaging). Because of our emphasis 
on particles with St n = 1, we also replace the upper limit of VJMR’s class 1 integral (k*) 
with the Kolmogorov scale k n . This simplification is, in fact, actually fairly good over the 
entire range of St L « 1, precisely because the contribution of eddies on smaller scales than 
k* (the class 3 eddies) is negligible. That is, the upper limit can be extended from k* to k n 
in general for mathematical simplicity without incurring significant error. Mathematically, 
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the upper limit could even be extended to infinity {eg., Volk et al. 1980), but the important 
role of the Reynolds number and of the Kolmogorov scale is then lost. Thus, 

V ?~ 2 £ m T^T. dk - «> 

Similarly, the generalized MMV expression for V* (their equation 6) can be simplified to 


K 2 



, 1 T tk/t s 



dk 


(5) 


for the particle size regime of interest here. As did VJMR, MMV note that the second 
integral of their equation (6) - the class 3 eddy contribution - is negligible for small particles, 
so we retain only the first integral of their equation (6). We again simplify the upper limit 
of integration in the remaining integral for the nominal St, } ~ 1 case where k* ~ k v >> ki. 
We validate this by comparing our results with those of MMV (section 1.7). 

The result for V 2 was plotted, but not stated explicitly, by VJMR and MMV (figure 1 in 
both papers), and explicitly derived by Cuzzi et al. (1993; Appendix B): V 2 = V 2 /(l + St L ). 
It is simple to see why V 2 « V 2 in the limit St L « 1 and certainly for St° « 1, since 
t s << tk in equations (4) or (5) for nearly all k and overwhelmingly all E{k). This limit is 
appropriate fpr chondrule-and-CAI-sized particles even in the presence of their small vertical 
settling velocity - they diffuse nearly as well as a gas molecule, and do not “settle to the 
midplane” in even a very weakly turbulent nebula (Dubrulle et al. 1995, Cuzzi et al. 1996). 
The implications are discussed in section 3. However, V 2 and V 2 are not exactly equal, 
resulting in a small, but very important, relative energy of motion V pg , giving the velocity 
with which particles move through the gas and encounter tiny (micron-sized) dust grains. 


1.4 Particle velocities relative to the gas 

The average relative velocity magnitude between a particle and the turbulent gas is V ps . 
VJMR make use only of the spatial frequency components of this quantity, which they refer to 
35 K-e/ {k) (their equation 15). Practically speaking, however, a particle will instantaneously 
sense all eddy contributions as one V pg ; we obtain this by merely integrating VJMR equation 
(15) over k. Considering only the part of the expression relevant for St v ~ 1 (that for 
k* > k L ), neglecting any systematic velocity, and again letting k* & k v » k L , the second 
integral vanishes and we obtain 



For this n — 0 case treated by VJMR, it can be easily verified using equations (4) and 
(6) that 

V p 2 9 + V p 2 = 2 f kv E(k)dk = V?. (7) 

J ki 

However, this useful result is true independent of n. It may also be obtained by Fourier 
transform solution of the forcing equations in temporal frequency ( u ) space, where the 
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energy spectrum of gas velocity fluctuations E g (u), particle velocity fluctuations E p ( u), and 
relative velocity fluctuations E pg (u) are related by 

E p (uj) = E g {u)/( 1 + t 2 s u 2 ) and E pg (ui) = t 2 s u 2 E p {uj). (8) 


This approach can be traced to Csanady (1963); it is also described by Hinze (1975, chapter 
5), Meek and Jones (1973), and Squires (1990, sections 4.2 and 4.5.1)). The E p solution 
was also derived in this way by Cuzzi et al. (1993, Appendix B). It is also clear then that 
E pg (uj ) + E p (u>) = Eg(u), essentially the same result as equation (7) above. Finally, we have 
directly verified equation (7) in our numerical simulations. 

Using this general relationship, we can extend the results of MMV to obtain V 2 g for their 
more generalized gas velocity autocorrelation functions (they only present results for V 2 ). 
That is, using equations (1) and (5), 


K = * - v; -2 £ Em^j-r'dk = 2 £ *<*> 


1 

1 + h/ts 


n+1 


dk. 


(9) 


We will use equations (5) and (9), with assumed inertial range expressions for E(k), to 
derive analytical expressions for V p and V pg of hypothetically “chondrule-like” (ie., St n ss 1) 
particles as functions of their size and the turbulent Reynolds number. 


1.5 Relative velocities between particles of similar sizes 


Expressions for V pp (VJMR Appendix C and equation 19; MMV equations 7 and 8) are more 
cumbersome, but respond nicely to certain simplifying assumptions. The full expression for 
V pp for two particles of equal size is (changing notation slightly from MMV equation 9, and 
allowing for a finite Kolmogorov scale): 




t s \ [ ^ , nt s h{x) 

tT+Tj g(x) + TT u 


dk, 


( 10 ) 


where g(x) = tan -I (x)/x and Mx) = 1/(1 + X 2 )’ The parameter x of VJMR and MMV is 
small in our regime of interest: 


K-ei(fc)ff(*F*:) _ Vrel{k)t s ~ V T e j(fc) 

ts + tk v(k)(t s + t k ) 2 v(k) 


( 11 ) 


since in the very limited range of k over which the integral is done, t s « t k . 1 In fact y << 1 
over most of the integral where t s « t k , so the functions g(x) and Hx) are ~ 1 or perhaps 
as small as a fraction of order unity; thus 


v 1 ~ 


rj 

/ m 

Jk> 


1 - 


t s N 


t s + tfc 


\n+n rk g 

j dk = 4 J E(k ) 


1 

1 + tk/t s 


n+l 


dk. (12) 


The integrand is identical to that for V 2 , but the integral has different limits which make it 
clear that only the eddies faster than t s can perturb identical particles into having incoherent 
relative velocities. 

Rn the above equation, the mathematical generalization of V pg by VJMR and MMV to its k-th com- 
ponents V ret {k) momentarily reappears. However, it is true in general, at any spatial frequency, that the 
particle-gas relative velocity is less than, or at most equal to, the gas velocity itself. 
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1.6 Scaling relations 

Recall that for the gas, 

U = m/v(k) = (L/V g )(k/k L )-W = ti(k/k L )-W (13) 

(Cuzzi et al. 2001). In equation (13) we have made the usual identification of V g with the 
largest scale eddy L. For the particles, 

St L = (k s /k L )~ 2/3 (14) 

and 

£ = (*/*«) 2/3 = ^(k/k L ) 2/3 = St L {k/k L ) 2/ 3 , (15) 

Note that if we restrict our attention to particles with St v = t s /t g « 1, then their Stokes 
number referred to the integral scale automatically becomes 

St L = t s /t L = t v /t L = (kr,/k L )~ 2/3 = (Re 3/4 )~ 2/3 = Re~ 1/2 . (16) 

The last substitution of (k^/ki) = Re 3 ^ 4 , where Re = LV g /v is the flow Reynolds number, 
with v being the molecular kinematic viscosity, is a direct consequence of the definitions of 
the Kolmogorov scale, the energy dissipation rate, and the Reynolds number (Tennekes and 
Lumley 1972). This relation can be obtained without any reference at all to the Kolmogorov 
spectrum but merely using scaling arguments relating to t L and t v . 2 Re is related to astro- 
physical “a”-models of the protoplanetary nebula by Re — acH/v with c = sound speed 
and H = nebula vertical scale height (Cuzzi et al. 2001). 


1.7 Final expressions for V pg and V pp 

Substituting the scaling relations from above for t s /t k , equation (9) for V pg becomes 



f -.. 1 

\1 + tk/t s 


n+l 


dk = 2 



St L 

St L + (k/k L )~ 2 / 3 


n+l 

dk. 


(17) 


We use the normalization (equation 1) to write E(k) = {V 2 /3k L )(k/k L )~ 5 / 3 , and change 
integration variable to x — k/k L , leaving 


V 2 

P9 


2V 2 rR* 3 '*/ Sh \ n+1 
3 Ji \StL + x~ 2 / 3 ) 


x-V 3 dx. 


(18) 


where in the upper limit we have substituted k p /ki = Re 3 ^ 4 from the scaling relations. 
Closed form solutions for equation (18) can be obtained for n = 0 or 1. For example, for 
n = 1 the result of the integral is 


St L 

1 

— V 2 

StHRe 1 / 2 - 1) 

.1 + Stix 2 / 3 . 

Re 3 / 4 

- V 9 

[(St L + l)(St L ReW + l) m 


(19) 


2 Let the energy dissipation rate be e. Then t — V 2 /ti = V g /L where the first expression defines 
and the last expression defines L. Also t v = (v/e) 1 / 2 and t? = (i^/e) 1 / 4 . Solving gives t L /t„ = Re 1 ! 2 and 
t i/L = Re 3 /*. 
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For n — 0 the result of the integral is: 



St L In 


(Re^il + StM 

{ Re l /*St L + 1 )_ 


(20) 


These results make it quite easy to predict both the magnitude and the St^ dependence of 
V pg for arbitrary nebula turbulent intensity. 

We solve equation (12) for V pp in a similar fashion to the solution for V pg above, to obtain 


V 2 = t Yl f k,,/kL ( 2gt L x~ 7 / 3 4 X- 9 / 3 \ 

pp 3 Jk{u)ik L \Stl + 2 St L x-W 4- x-^j 


( 21 ) 


As before, the upper integration limit is k v /k L = J?e 3/4 . For the lower limit, k* jk L — 
k(t s )/k L = ( t s /t L )~ 3 / 2 = St~ L 3/2 from the scaling relations. The closed form analytic solution 
of this integral is: 


V 2 = 2V 2 

PP 9 


X 


- 2/3 


[l + St L xV\ 


St~ 3/2 


= 2V‘ 


r st, 


i 


Re 3 /' * 


L 2 


St L Re 4 Re l > 2 


( 22 ) 


The n = 0 form of the solution is somewhat less useful, and we note it without expanding it 
as it will not be used further. 


V 2 

pp 


2V„ 


St i In 


I 4 St L x 2/ ' 3 ' 

^. 2/3 


n Re 3 / 4 


x 2/3 


Si 


- 3/2 


1.8 Detailed comparisons with the models of Markiewicz et al. 

In addition to developing the analytical expressions discussed and applied in the paper, 
we also developed a detailed numerical model following the prescriptions of MMV exactly 
(but with a generalized turbulent energy spectrum). This was needed both to evaluate 
their theoretical approach in the context of our numerical simulations of turbulence (section 
2), which have a non-Kolmogorov spectrum and low Reynolds number compared to nebula 
applications, and to assess the validity of our analytical approximations. The numerical 
model of MMV is no longer in active use (W. Markiewicz, personal communication 2002), 
so we digitized their V pp results (their figure 5) to facilitate comparisons. As seen in figure 
4, our full numerical model for V pp (solid curves) agrees very well with their results for V pp 
(long dashed curves). In figure 2 we also show our results for V pg , not presented by VJMR 
or MMV, as obtained by integrating MMV equation 4 over all spatial frequencies. Note that 
we, and MMV, both use the appropriate form of R(t,t',k ) (ie., that for the correct choice 
of n; section 1.2) for these calculations. 

The most striking feature of the results, first noted by MMV, is that V pp very quickly falls 
to zero for particles with St n < 1 (i.e. St L < Re~ l/2 , as shown in the scaling relations of sec- 
tion 1.6 above) because there is no more energy in faster eddies to provide relative velocities 
to such particles. This does not happen to V pg , because eddies on all scales contribute. Also 
note that V p and V pp decrease for large particles (Sti > 1), as fewer eddies can effectively 
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st L = t./t L 


Figure 2: Comparison of our numerical version of the full MMV model for n = 0 (light 
curves) and n = 1 (heavy curves), along with digitized results for V pp from MMV (dashed 
curves, their figure 5, for n = 1). Three different Re are shown: (a) 10 4 , (b) 10 7 , and (c) 
10 9 . The dash-dot curves are for V p , which has the same shape for all three Re. V vg is shown 
in the two sets of dotted curves and V pp in the two sets of solid curves. Note that the n = 0 
values of V pg (light dotted curves) are considerably (3-4 times) higher than the preferred 
n — 1 values (heavy dotted curves), and the SjT-dependence of V pg , for n = 0, never gets 
much above 0.5, whereas for n = 1 a linear dependence is seen for St L < Re~ 1/2 . As in 
figures 4 and 5, vertical hash marks correspond to St L = Re~ l/2 for the three values of Re. 
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couple to particles with such long stopping times. Naturally, V pg simply approaches V g for 
these large particles. 

Upon comparing our original analytical results (equations 19 and 22) with our full nu- 
merical model and the MMV results, we found some small quantitative discrepancies at the 
order unity level, as might be expected. The responsible approximations were easily identi- 
fied. First, we approximated the boundary between class 1 and class 3 eddies by t s = t(k*) 
rather than the more complete equation 9 of VJMR and equation 4 of MMV, which ob- 
tains the relevant eddy frequency in the moving frame of the particle and involves V re i(k). 
Comparison of the two criteria revealed that, to a very good approximation, the criterion 
t s = t(k*) gives a value of k* that is too large by a factor close to 2 (figure 3). So, after 
this “calibration”, we merely decrease the lower limit of integration in our equation (22) by 
a factor of 2. Second, even after this correction, our values of V vp are about 20% high. This 
is easily ascribed to our approximation that g(x) and h(x) are equal to unity throughout 
the entire range of A;; in fact, they are tens of percent smaller than unity over some part 
of this range, depending on the value of St L . Empirically, this is corrected by multiplying 
our analytical expression for V pp by a constant factor of 0.8. With these two simple adjust- 
ments, each correcting a known oversimplification, our analytical expression for V pp achieves 
very good good agreement with the MMV results, and with our own full numerical model, 
over the relevant range of St L < 0.1 or so. There appears to be no reason to make such 
refinements to our analytical expression for V pg (equati&n 19), because our approximations 
are better justified and the agreement with MMV acceptable. 


1.9 Numerical refinements to the model 

With insights gained from comparison of our numerical and analytical models, we have made 
two small adjustments to equation (22) for V pp which correct for two of our approximations. 
Equation (22) is multiplied by a factor of 0.8, and the upper integration limit (Sfjj 3 ^) is 
divided by two, so the first term in the final expression changes from St L / 2 to St L j 1.03 « 
St L . The approximations entering into our expression for V pg are better, so no correction is 
applied. The final equation for V pp is then 


V 2 = 1.6 V 2 

PP S 


St l — 


1 

St^Re ~+~Re^\ 


(23) 


The results of equations (19) and (23) (the preferred and adjusted n — 1 forms), normalized 
by V 9 , are shown in figure 4 for the same three values of Re as in MMV, and in closeup 
form in figure 5. 

As shown by MMV (their figure 2), and as seen previously in our figure 2, the falloff 
of Vpp is extremely steep for St^ < 1 (i.e. St L < Re~ 1/2 as shown in the scaling relations 
of section 1.6 above) because there is no more energy in faster eddies to provide relative 
velocities to such particles. 


1.10 Simplification of analytically determined velocity expressions: 

Equations (19) and (23) - for the preferred n = 1 case - are readily simplified in different 
limits of interest. It is simply shown by retaining leading terms that equation (19) for V pg 
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Stt 


Figure 3: Correction of our approximation k* & St by a factor of two (dash-dot line) 
which brings it into excellent agreement (in our range of validity St L < 0.1) with the ex- 
act numerical solution for k * , shown for Re = 10 4 , 10 7 , and 10 9 , computed using the full 
VJMR/MMV expression. Only very close to St^ = 1 does our approximation deviate slightly; 
notice the tiny tail at St L = 6 x 10~ 3 , k* = 10 3 , which is the Kolmogorov scale for Re = 10 4 . 
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Figure 4: V pg (St L ) (dotted; equation 19) and V pp (St L ) (solid; equation 23) for Re = 10 4 (a), 
10 7 (b), and 10 9 (c). The digitized results of MMV (theif figure 5) for V pp , for the same three 
values of Re, are shown by the dashed lines. Our V pp expression is invalid for St L > 0.1 or 
so (see text). 



Figure 5: A closeup plot of V pg (dotted), V pp (solid), and the digitized MMV results (their 
figure 5) for V pp (dashed, see Appendix) all for the Re = 10 7 case. The dash-dot line has 
slope 1/2. The vertical short dashed line indicates St ^ = 1, where St L = Re~ 1/2 \ here, 
V pg oc St° L 75 . 
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results in three separate regimes: V P 9 ~ v g for St L > 1, V pg oc St l L /2 for Re 1/2 < St L « 1, 
and V pg oc St L Re 1 ^ for Sti < Re ~ 1 - /2 . This is confirmed by inspection of figures 2 and 4. In 
the special case of St v = 1, or Sti = Re ~ 1/l2 , equation (19) reduces directly to 

v,,{stn " = v ’~yT = m ' 14 {'£h) ’ (24) - 

where we have substituted V g = ca 1/2 (Cuzzi et al. 2001). This Re-dependence, which also 
applies for St v < 1 in general, quite naturally explains a result we obtained empirically from 
our numerical models over a range of Re much smaller than nebula values, namely that 
Vpg/Vg K R e (Cuzzi et al. 1998). By contrast, it is similarly shown from equation (20) 

that the SY L -dependence 0 f V pg for the older n = 0 case continues the St 1 / 2 dependence to 
arbitrarily small St L . 

These results are also consistent with arguments in Cuzzi et al. (1993, Appendix B; 
A. Dobrovolskis, personal communication). Expand and time-average the instantaneous 
quantity < (V p - V g ) 2 > to obtain < V pg V pg >=< V p V p > + < V g V g > -2 < V p V g >. 
Substituting from Cuzzi et al. (1993, equation Bll) we find < V p V p > = < V p V g >=< V g V g > 
/(1 + Sti), leading to V pg = (Sti/( 1 + Sti)) 1 ^ 2 V g , which reaches the same limits as equation 
(19) except for particles with t s < t v , or St v < 1, because the integral in its derivation 
(equation Bll of Cuzzi et al. 1993) extends to infinite eddy frequency. 

Thus, unless t s < t v ( Sti < Re -1 / 2 ), the particle-gas relative velocity in turbulence is 
generally proportional to y/Sti for small Sti- The steeper dependence of V pg on Sti and 
St v is restricted (in turbulence) to particles with St n < 1. That is, evidence for a more 
nearly linear dependence of V pg on r, if the environment was turbulent, would imply that 
the particles in question were St^ < 1 particles. This new result derives directly from the 
use of the n = 1 gas velocity autocorrelation function. The primary qualitative change is in 
the particle size dependence of V pg for particles with St v < 1. We address the significance 
of this in more detail in a forthcoming paper (Cuzzi 2002b). 

Finally, using equation (23) for V pp , we get 

Vpp{St n = 1) = VoSVgRe^ 4 = 1.2QV pg , (25) 

where we used equation (24) for V pg . 

2 Comparison with numerical results 

In this section we compare numerical results from our full 3D Lagrangian particle-gas model 
(Hogan et al. 1999) with full numerical calculations using our implementation of MMV 
(sections 1.8- 1.9). We present particle velocities relative to the computational box (V p ), and 
relative to the local fluid velocity (V^), as obtained from our simulations. These velocities 
are defined as RMS spatial averages over all particles in a single snapshot, or V = (< 04 - < 

Vx >) > + < ( Vy — < Vy >) 2 > + < (V z — < V z >) 2 >) 1/,2 ) where V represents V p or V pg at 
the location of each particle, and <> is the averaging operator < ... >= where 

N p is the number of particles in a single snapshot. Of course, < V > is very close to zero 
for both these quantities since there is no mean flow in our simulations. 
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st L 

Figure 6: V p vs. St L obtained from our direct simulations compared with MMV predictions 
for models n = 0 (solid line) and n = 1 (dashed line). AH velocities have been normalized by 
the RMS fluid velocity V g . Results are shown for three different Re] the St L values for each 
point are defined relative to a large eddy time based on energy dissipation 2 , which varies 
with Re for our numerical calculations. When they are defined relative to a constant large 
eddy time, as are our analytical models and the MMV models shown in figure 2, points 
and models for all Re collapse onto the same curve as seen in figure 2. The n = 1 MMV 
prediction is clearly a better fit to the numerically simulated velocities, regardless of the 
choice of normalization timescale. 

This spatial averaging approach is equivalent to the temporal averaging implicit in the 
MMV model, because of the ergodic principle that equates temporal and spatial averaging 
under suitable conditions. In our case, the conditions are satisified because our integral 
length scale L is small compared to the spatial period of the computational domain, for all 
Re. 

The case of V pp is more complicated, as the results depend on the proximity region chosen 
for “neighboring” particles. For the most useful comparisons with the predictions of MMV 
and VJMR, and with the expected uses of this quantity in mind, the region over which 
particle neighbors are selected should be as small as possible - less than r] certainly - and 
here we run into sampling errors. Perhaps most important, the deviation of our model energy 
spectrum from a Kolmogorov spectrum is significant (eg. Squires and Eaton 1990), and V pp is 
much more sensitive to the details of the high-spatial-frequency end of the energy spectrum 
than either V p or V pg . Since the main purpose of these calculations is to verify numerically 
the preference for the n = 1 autocorrelation function in an independent way from the direct 
comparison shown in figure 1, and because this case is already well made by the V p and V pg 
plots, we present no comparisons for V pp . 

Figures 6-9 show that the n — 1 autocorrelation function provides a much better 


15 





st L 


Figure 8: V pg obtained from our direct simulations, compared with predictions of the MMV 
models with n = 0 (solid line) and n = 1 (dashed line) vs St j J . All curves have been 
normalized by the RMS fluid velocity V g 





Figure 9: The same data shown in figure 8, but plotted against St v . 

fit to both V p and V pg than the n = 0 version. For V pg , the fits of the MMV theory to 
our simulations are less perfect than for V p . We can see several possible explanations for 
this. For instance, the mathematically simple form adopted for the n = 1 autocorrelation 
function is not a perfect fit to the actual numerically determined one (figure 1), by about the 
correct fractional amount. Also, we have emphasized that the correct velocity autocorrelation 
function to use is that along a particle trajectory (Meek and Jones 1973), and this function 
is actually somewhat size dependent even over the range Stp ~ 1 (see, eg., Squires 1990, 
figure 4-23). Finally, because of the deviation of our turbulent kinetic energy spectrum from 
an inertial range, some of the definitions of eddy times used in the MMV theory might 
be inappropriate. It would not be surprising for V pg to be more sensitive to these small 
deviations than V v (compare figures 6 and 8, or 7 and 9). In spite of the small deviations in 
V vg , the combination of the direct comparisons of the autocorrelation functions themselves 
(figure 1), and the comparison of the velocities derived using them (figures 6-9) makes it 
clear that the n = 1 autocorrelation function is the best choice. 


3 Summary and conclusions 

We present theoretical and numerical results which describe the turbulence-driven velocities 
of particles in the St L « 1 size regime which might characterize chondrules and similar 
sized particles. We numerically verify the general approach of VJMR as modified by MMV, 
and verify in two different ways the intuitive preference of MMV for an n = 1 gas velocity 
autocorrelation function - at least along the trajectories of St v « 1 particles. We find 
theoretically that the n = 1 autocorrelation function leads to a particle-gas relative velocity 
function that approaches linear dependence on particle size for particles in the Stg « 1 
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regime, and becomes and remains linear for arbitrarily small sizes. This is quite a different 
result than predicted by the original VJMR n — 0 expressions. We derive simple analytic 
expressions for V p , V pg , and V pp (the latter, for comparable size particles only) for arbitrary 
levels of nebula intensity, as characterized by its Reynolds number Re or its corresponding 
“a”. In a separate paper (Cuzzi 2002b) we will present some implications of these results 
for meteoritics. 
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Table I. List of symbols 


Parameter 

c 

C 

E(k) 

H 

k 

k L 

k rj 

L 

r 

R 

Re 

St i 

Str, 

ts 

tk 

tL 

tr) 

V 9 

v p 

V P9 

v pp 

a 

t 

V 

v 

Vj' 

u 

Pg 

Ps 


De finition 

gas molecule thermal speed 

particle concentration factor 

turbulent kinetic energy at wavenumber k 

nebula vertical scale height 

eddy wavenumber 

wavenumber of largest eddy 

wavenumber of Kolmogorov scale eddy 

integral or largest scale in turbulent energy spectrum 

particle radius 

gas velocity autocorrelation function 

flow Reynolds number 

Stokes number relative to largest eddy 

Stokes number relative to Kolmogorov scale eddy 

stopping time of particle due to gas drag 

overturn time of eddy with wavenumber k 

overturn time of largest eddy 

overturn time of Kolmogorov scale eddy 

gas turbulent velocity (large eddy) 

particle random velocity in inertial space 

relative velocity between particles and gas 

relative velocity between particles B 

nebula viscosity parameter; Re — acH/ v 

dissipation of turbulent kinetic energy 

Kolmogorov scale 

molecular kinematic viscosity 

turbulent kinematic viscosity 

eddy temporal frequency 

gas mass density 

solid particle mass density 
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